Gerrymandering Notebook

library(sf)          # classes and functions for vector data

Linking to GEOS 3.8.0, GDAL 3.0.4, PROJ 6.3.1; sf_use_s2() is TRUE

library(spData)        # load geographic data

To access larger datasets in this package, install the spDataLarge package with: `install.packages('spDataLarge', repos='https://nowosad.github.io/drat/', type='source')`

library(dplyr)

Attaching package: ‘dplyr’

The following objects are masked from ‘package:stats’:

    filter, lag

The following objects are masked from ‘package:base’:

    intersect, setdiff, setequal, union

library(tmap)    # for static and interactive maps

Registered S3 method overwritten by 'htmlwidgets':
  method           from         
  print.htmlwidget tools:rstudio

library(ggplot2) 
library(ggsci)
library(tools)
library(readr)

library(foreach)
library(doParallel)

Loading required package: iterators
Loading required package: parallel

registerDoParallel(detectCores())
getDoParWorkers()

[1] 96

load("data/ca.data.RData")
dem=tm_shape(ca.data) +tm_fill(col = "dist",palette = "viridis",legend.show=FALSE,alpha=0.4)+
  tm_dots(size = "party_dem",col="blue",alpha=0.5,legend.show=FALSE)+
  tm_layout(title = "California", frame = FALSE, bg.color = NA)
rep=tm_shape(ca.data) +tm_fill(col = "dist",palette = "viridis",legend.show=FALSE,alpha=0.4)+
  tm_dots(size = "party_rep",col="red",alpha=0.5,legend.show=FALSE)+
  tm_layout(title = "California", frame = FALSE, bg.color = NA)
black=tm_shape(ca.data) +tm_fill(col = "dist",palette = "viridis",legend.show=FALSE,alpha=0.4)+
  tm_dots(size = "eth1_aa",col="black",alpha=0.5,legend.show=FALSE)+
  tm_layout(title = "California", frame = FALSE, bg.color = NA)
white=tm_shape(ca.data) +tm_fill(col = "dist",palette = "viridis",legend.show=FALSE,alpha=0.4)+
  tm_dots(size = "eth1_eur",col="yellow",alpha=0.5,legend.show=FALSE)+
  tm_layout(title = "California", frame = FALSE, bg.color = NA)
hisp=tm_shape(ca.data) +tm_fill(col = "dist",palette = "viridis",legend.show=FALSE,alpha=0.4)+tm_dots(size = "eth1_hisp",col="brown",alpha=0.5,legend.show=FALSE)+tm_layout(title = "California", frame = FALSE, bg.color = NA)
asian=tm_shape(ca.data) +tm_fill(col = "dist",palette = "viridis",legend.show=FALSE,alpha=0.4)+tm_dots(size = "eth1_esa",col="orange",alpha=0.5,legend.show=FALSE)+tm_layout(title = "California", frame = FALSE, bg.color = NA)


tmap_arrange(black,white,hisp,asian,dem,rep,nrow = 2,ncol=4)

load("data/tx.data.RData")
dem=tm_shape(tx.data) +tm_fill(col = "dist",palette = "viridis",legend.show=FALSE,alpha=0.4)+
  tm_dots(size = "party_dem",col="blue",alpha=0.5,legend.show=FALSE)+
  tm_layout(title = "Texas", frame = FALSE, bg.color = NA)
rep=tm_shape(tx.data) +tm_fill(col = "dist",palette = "viridis",legend.show=FALSE,alpha=0.4)+
  tm_dots(size = "party_rep",col="red",alpha=0.5,legend.show=FALSE)+
  tm_layout(title = "Texas", frame = FALSE, bg.color = NA)
black=tm_shape(tx.data) +tm_fill(col = "dist",palette = "viridis",legend.show=FALSE,alpha=0.4)+
  tm_dots(size = "eth1_aa",col="black",alpha=0.5,legend.show=FALSE)+
  tm_layout(title = "Texas", frame = FALSE, bg.color = NA)
white=tm_shape(tx.data) +tm_fill(col = "dist",palette = "viridis",legend.show=FALSE,alpha=0.4)+
  tm_dots(size = "eth1_eur",col="yellow",alpha=0.5,legend.show=FALSE)+
  tm_layout(title = "Texas", frame = FALSE, bg.color = NA)
hisp=tm_shape(tx.data) +tm_fill(col = "dist",palette = "viridis",legend.show=FALSE,alpha=0.4)+tm_dots(size = "eth1_hisp",col="brown",alpha=0.5,legend.show=FALSE)+tm_layout(title = "Texas", frame = FALSE, bg.color = NA)
asian=tm_shape(tx.data) +tm_fill(col = "dist",palette = "viridis",legend.show=FALSE,alpha=0.4)+tm_dots(size = "eth1_esa",col="orange",alpha=0.5,legend.show=FALSE)+tm_layout(title = "Texas", frame = FALSE, bg.color = NA)


tmap_arrange(black,white,hisp,asian,dem,rep,nrow = 2,ncol=4)

Legend labels were too wide. Therefore, legend.text.size has been set to 0.59. Increase legend.width (argument of tm_layout) to make the legend wider and therefore the labels larger.
Legend labels were too wide. Therefore, legend.text.size has been set to 0.59. Increase legend.width (argument of tm_layout) to make the legend wider and therefore the labels larger.

NA

library(ggplot2)
library(dplyr)
library(ggsci)
library(qte)

Registered S3 method overwritten by 'data.table':
  method           from
  print.data.table     

library(CVXR)

Attaching package: ‘CVXR’

The following object is masked from ‘package:dplyr’:

    id

The following object is masked from ‘package:stats’:

    power

load("data/power.data.RData")

###### My own Quantile treatment effect model


quantile=function(taus, y, wt){
 
  # minimize sum of check functions
  
  q <- Variable(1)
  quant_loss <- function(u, tau) { 0.5 * abs(u) + (tau - 0.5) * u }
  solutions <- sapply(taus, function(tau) {
    obj <- sum(wt*quant_loss(y-q,tau))
    prob <- Problem(Minimize(obj))
    ## THE OSQP solver returns an error for tau = 0.5
    solve(prob, solver = "ECOS")$getValue(q)
  })
  return(solutions)
}

QTE=function(dataset, taus){
  
  
  #pscore.reg <- glm(T ~ , data=dataset, family=binomial)
  #pscore <- fitted(pscore.reg)
  
  q.trt=quantile(taus=taus, y=dataset$Y, wt=dataset$T*dataset$eth1_hisp)
  q.control=quantile(taus=taus, y=dataset$Y, wt=(1-dataset$T)*dataset$eth1_hisp)
  return(q.trt-q.control)
  
}

#unique(power.data$STATE)
power.data%>%
  mutate(black.power=black.power*100,
         white.power=white.power*100,
         hisp.power=hisp.power*100,
         asian.power=asian.power*100)%>%
  mutate(A=ifelse(STATE %in% c("Arizona","California","Hawaii","Idaho","Montana","New Jersey","Washington","Alaska","Arkansas","Colorado","Missouri","Ohio","Pennsylvania","Iowa"),
                  "Nonpartisan or bipartisan commissions",
                  ifelse(STATE %in% c("Delaware","Illinois","Maryland","Massachusetts","Nevada","New Mexico","New York","Oregon","Rhode Island"),
                         "Democratic",
                         ifelse(STATE %in% c("Alabama","Florida","Georgia","Indiana","Kansas","Kentucky","Mississippi","New Hampshire","North Carolina","North Dakota","Oklahoma","South Carolina","South Dakota","Tennessee","Texas","Utah","West Virginia","Wyoming"),"Republican",""))))%>%filter(A !="")%>%
  ggplot()+
  geom_density(aes(x=white.power,y=..scaled..,weight=eth1_eur,group=A,col=A))+
  scale_color_startrek()+
  ggtitle("White voter power")

power.data%>%
  mutate(black.power=black.power*100,
         white.power=white.power*100,
         hisp.power=hisp.power*100,
         asian.power=asian.power*100)%>%
  mutate(A=ifelse(STATE %in% c("Arizona","California","Hawaii","Idaho","Montana","New Jersey","Washington","Alaska","Arkansas","Colorado","Missouri","Ohio","Pennsylvania","Iowa"),
                  "Nonpartisan or bipartisan commissions",
                  ifelse(STATE %in% c("Delaware","Illinois","Maryland","Massachusetts","Nevada","New Mexico","New York","Oregon","Rhode Island"),
                         "Democratic",
                         ifelse(STATE %in% c("Alabama","Florida","Georgia","Indiana","Kansas","Kentucky","Mississippi","New Hampshire","North Carolina","North Dakota","Oklahoma","South Carolina","South Dakota","Tennessee","Texas","Utah","West Virginia","Wyoming"),"Republican",""))))%>%filter(A !="")%>%
  ggplot()+
  geom_density(alpha=0.1,aes(x=black.power,y=..scaled..,weight=eth1_aa,group=A,col=A))+scale_color_startrek()+
  ggtitle("Black voter power")

 
power.data%>%
  mutate(black.power=black.power*100,
         white.power=white.power*100,
         hisp.power=hisp.power*100,
         asian.power=asian.power*100)%>%
  mutate(A=ifelse(STATE %in% c("Arizona","California","Hawaii","Idaho","Montana","New Jersey","Washington","Alaska","Arkansas","Colorado","Missouri","Ohio","Pennsylvania","Iowa"),
                  "Nonpartisan or bipartisan commissions",
                  ifelse(STATE %in% c("Delaware","Illinois","Maryland","Massachusetts","Nevada","New Mexico","New York","Oregon","Rhode Island"),
                         "Democratic",
                         ifelse(STATE %in% c("Alabama","Florida","Georgia","Indiana","Kansas","Kentucky","Mississippi","New Hampshire","North Carolina","North Dakota","Oklahoma","South Carolina","South Dakota","Tennessee","Texas","Utah","West Virginia","Wyoming"),"Republican",""))))%>%filter(A !="")%>%
  ggplot()+
  geom_density(alpha=0.1,aes(x=hisp.power,y=..scaled..,weight=eth1_hisp,group=A,col=A))+scale_color_startrek()+
  ggtitle("Hisp voter power")

power.data%>%
  mutate(black.power=black.power*100,
         white.power=white.power*100,
         hisp.power=hisp.power*100,
         asian.power=asian.power*100)%>%
  mutate(A=ifelse(STATE %in% c("Arizona","California","Hawaii","Idaho","Montana","New Jersey","Washington","Alaska","Arkansas","Colorado","Missouri","Ohio","Pennsylvania","Iowa"),
                  "Nonpartisan or bipartisan commissions",
                  ifelse(STATE %in% c("Delaware","Illinois","Maryland","Massachusetts","Nevada","New Mexico","New York","Oregon","Rhode Island"),
                         "Democratic",
                         ifelse(STATE %in% c("Alabama","Florida","Georgia","Indiana","Kansas","Kentucky","Mississippi","New Hampshire","North Carolina","North Dakota","Oklahoma","South Carolina","South Dakota","Tennessee","Texas","Utah","West Virginia","Wyoming"),"Republican",""))))%>%filter(A !="")%>%
  ggplot()+
  geom_density(alpha=0.1,aes(x=asian.power,y=..scaled..,weight=eth1_esa,group=A,col=A))+scale_color_startrek()+
  ggtitle("Asian voter power")

 
#table(power.data$STATE,power.data$A)

#unique(power.data$STATE)
power.data=power.data%>%
  mutate(black.power=black.power*100,
         white.power=white.power*100,
         hisp.power=hisp.power*100,
         asian.power=asian.power*100)%>%
  mutate(A=ifelse(STATE %in% c("Arizona","California","Hawaii","Idaho","Montana","New Jersey","Washington","Alaska","Arkansas","Colorado","Missouri","Ohio","Pennsylvania","Iowa"),
                  "Nonpartisan or bipartisan commissions",
                  ifelse(STATE %in% c("Delaware","Illinois","Maryland","Massachusetts","Nevada","New Mexico","New York","Oregon","Rhode Island"),
                         "Democratic",
                         ifelse(STATE %in% c("Alabama","Florida","Georgia","Indiana","Kansas","Kentucky","Mississippi","New Hampshire","North Carolina","North Dakota","Oklahoma","South Carolina","South Dakota","Tennessee","Texas","Utah","West Virginia","Wyoming"),"Republican",""))))%>%
  filter(A %in% c("Republican","Nonpartisan or bipartisan commissions"))%>%mutate(T=ifelse(A=="Republican",1,0))
#table(power.data$STATE,power.data$A)

taus=c(0.05,0.1, 0.25, 0.5, 0.75, 0.90, 0.95)
q.trt=quantile(taus=taus, y=power.data$black.power, wt=power.data$T*power.data$eth1_aa)

Warning in size(const) : NAs introduced by coercion to integer range
Warning in size(const) : NAs introduced by coercion to integer range
Warning in size(const) : NAs introduced by coercion to integer range
Warning in size(const) : NAs introduced by coercion to integer range
Warning in size(const) : NAs introduced by coercion to integer range
Warning in size(const) : NAs introduced by coercion to integer range
Warning in size(const) : NAs introduced by coercion to integer range
Warning in size(const) : NAs introduced by coercion to integer range
Warning in size(const) : NAs introduced by coercion to integer range
Warning in size(const) : NAs introduced by coercion to integer range
Warning in size(const) : NAs introduced by coercion to integer range
Warning in size(const) : NAs introduced by coercion to integer range
Warning in size(const) : NAs introduced by coercion to integer range
Warning in size(const) : NAs introduced by coercion to integer range

q.control=quantile(taus=taus, y=power.data$black.power, wt=(1-power.data$T)*power.data$eth1_aa)

Warning in size(const) : NAs introduced by coercion to integer range
Warning in size(const) : NAs introduced by coercion to integer range
Warning in size(const) : NAs introduced by coercion to integer range
Warning in size(const) : NAs introduced by coercion to integer range
Warning in size(const) : NAs introduced by coercion to integer range
Warning in size(const) : NAs introduced by coercion to integer range
Warning in size(const) : NAs introduced by coercion to integer range
Warning in size(const) : NAs introduced by coercion to integer range
Warning in size(const) : NAs introduced by coercion to integer range
Warning in size(const) : NAs introduced by coercion to integer range
Warning in size(const) : NAs introduced by coercion to integer range
Warning in size(const) : NAs introduced by coercion to integer range
Warning in size(const) : NAs introduced by coercion to integer range
Warning in size(const) : NAs introduced by coercion to integer range

black.qte=q.trt-q.control
#plot(c(0.05,0.1, 0.25, 0.5, 0.75, 0.90, 0.95),black.qte,type="b")

q.trt=quantile(taus=taus, y=power.data$white.power, wt=power.data$T*power.data$eth1_eur)

Warning in size(const) : NAs introduced by coercion to integer range
Warning in size(const) : NAs introduced by coercion to integer range
Warning in size(const) : NAs introduced by coercion to integer range
Warning in size(const) : NAs introduced by coercion to integer range
Warning in size(const) : NAs introduced by coercion to integer range
Warning in size(const) : NAs introduced by coercion to integer range
Warning in size(const) : NAs introduced by coercion to integer range
Warning in size(const) : NAs introduced by coercion to integer range
Warning in size(const) : NAs introduced by coercion to integer range
Warning in size(const) : NAs introduced by coercion to integer range
Warning in size(const) : NAs introduced by coercion to integer range
Warning in size(const) : NAs introduced by coercion to integer range
Warning in size(const) : NAs introduced by coercion to integer range
Warning in size(const) : NAs introduced by coercion to integer range

q.control=quantile(taus=taus, y=power.data$white.power, wt=(1-power.data$T)*power.data$eth1_eur)

Warning in size(const) : NAs introduced by coercion to integer range
Warning in size(const) : NAs introduced by coercion to integer range
Warning in size(const) : NAs introduced by coercion to integer range
Warning in size(const) : NAs introduced by coercion to integer range
Warning in size(const) : NAs introduced by coercion to integer range
Warning in size(const) : NAs introduced by coercion to integer range
Warning in size(const) : NAs introduced by coercion to integer range
Warning in size(const) : NAs introduced by coercion to integer range
Warning in size(const) : NAs introduced by coercion to integer range
Warning in size(const) : NAs introduced by coercion to integer range
Warning in size(const) : NAs introduced by coercion to integer range
Warning in size(const) : NAs introduced by coercion to integer range
Warning in size(const) : NAs introduced by coercion to integer range
Warning in size(const) : NAs introduced by coercion to integer range

white.qte=q.trt-q.control
#plot(c(0.05,0.1, 0.25, 0.5, 0.75, 0.90, 0.95),white.qte,type="b")

q.trt=quantile(taus=taus, y=power.data$hisp.power, wt=power.data$T*power.data$eth1_hisp)

Warning in size(const) : NAs introduced by coercion to integer range
Warning in size(const) : NAs introduced by coercion to integer range
Warning in size(const) : NAs introduced by coercion to integer range
Warning in size(const) : NAs introduced by coercion to integer range
Warning in size(const) : NAs introduced by coercion to integer range
Warning in size(const) : NAs introduced by coercion to integer range
Warning in size(const) : NAs introduced by coercion to integer range
Warning in size(const) : NAs introduced by coercion to integer range
Warning in size(const) : NAs introduced by coercion to integer range
Warning in size(const) : NAs introduced by coercion to integer range
Warning in size(const) : NAs introduced by coercion to integer range
Warning in size(const) : NAs introduced by coercion to integer range
Warning in size(const) : NAs introduced by coercion to integer range
Warning in size(const) : NAs introduced by coercion to integer range

q.control=quantile(taus=taus, y=power.data$hisp.power, wt=(1-power.data$T)*power.data$eth1_hisp)

Warning in size(const) : NAs introduced by coercion to integer range
Warning in size(const) : NAs introduced by coercion to integer range
Warning in size(const) : NAs introduced by coercion to integer range
Warning in size(const) : NAs introduced by coercion to integer range
Warning in size(const) : NAs introduced by coercion to integer range
Warning in size(const) : NAs introduced by coercion to integer range
Warning in size(const) : NAs introduced by coercion to integer range
Warning in size(const) : NAs introduced by coercion to integer range
Warning in size(const) : NAs introduced by coercion to integer range
Warning in size(const) : NAs introduced by coercion to integer range
Warning in size(const) : NAs introduced by coercion to integer range
Warning in size(const) : NAs introduced by coercion to integer range
Warning in size(const) : NAs introduced by coercion to integer range
Warning in size(const) : NAs introduced by coercion to integer range

hisp.qte=q.trt-q.control
#plot(c(0.05,0.1, 0.25, 0.5, 0.75, 0.90, 0.95),hisp.qte,type="b")


q.trt=quantile(taus=taus, y=power.data$asian.power, wt=power.data$T*power.data$eth1_esa)

Warning in size(const) : NAs introduced by coercion to integer range
Warning in size(const) : NAs introduced by coercion to integer range
Warning in size(const) : NAs introduced by coercion to integer range
Warning in size(const) : NAs introduced by coercion to integer range
Warning in size(const) : NAs introduced by coercion to integer range
Warning in size(const) : NAs introduced by coercion to integer range
Warning in size(const) : NAs introduced by coercion to integer range
Warning in size(const) : NAs introduced by coercion to integer range
Warning in size(const) : NAs introduced by coercion to integer range
Warning in size(const) : NAs introduced by coercion to integer range
Warning in size(const) : NAs introduced by coercion to integer range
Warning in size(const) : NAs introduced by coercion to integer range
Warning in size(const) : NAs introduced by coercion to integer range
Warning in size(const) : NAs introduced by coercion to integer range

q.control=quantile(taus=taus, y=power.data$asian.power, wt=(1-power.data$T)*power.data$eth1_esa)

Warning in size(const) : NAs introduced by coercion to integer range
Warning in size(const) : NAs introduced by coercion to integer range
Warning in size(const) : NAs introduced by coercion to integer range
Warning in size(const) : NAs introduced by coercion to integer range
Warning in size(const) : NAs introduced by coercion to integer range
Warning in size(const) : NAs introduced by coercion to integer range
Warning in size(const) : NAs introduced by coercion to integer range
Warning in size(const) : NAs introduced by coercion to integer range
Warning in size(const) : NAs introduced by coercion to integer range
Warning in size(const) : NAs introduced by coercion to integer range
Warning in size(const) : NAs introduced by coercion to integer range
Warning in size(const) : NAs introduced by coercion to integer range
Warning in size(const) : NAs introduced by coercion to integer range
Warning in size(const) : NAs introduced by coercion to integer range

asian.qte=q.trt-q.control
#plot(c(0.05,0.1, 0.25, 0.5, 0.75, 0.90, 0.95),asian.qte,type="b")

data.frame(q=rep(c(0.05,0.1, 0.25, 0.5, 0.75, 0.90, 0.95),4),qte=c(black.qte,white.qte,hisp.qte,asian.qte),race=rep(c("Black","White","Hispanic","Asian"),each=7))%>%
  ggplot(aes(x=q,y=qte,col=race))+
  geom_point()+
  geom_line()+
  scale_color_startrek()